Metallic nanograins: spatially nonuniform pairing induced by quantum confinement 
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It is well-known that the formation of discrete electron levels strongly influences the pairing in 
metallic nanograins. Here we focus on another effect of quantum confinement in superconduct- 
ing grains that was not studied previously, i.e., spatially nonuniform pairing. This effect is very 
significant when single-electron levels form bunches and/or a kind of shell structure: in highly 
symmetric grains the order parameter can exhibit variations with position by an order of magni- 
tude. Nonuniform pairing is closely related to a quantum-confinement induced modification of the 
pairing-interaction matrix elements and size-dependent pinning of the chemical potential to groups 
of degenerate or nearly degenerate levels. For illustration we consider spherical metallic nanograins. 
We show that the relevant matrix elements are as a rule enhanced in the presence of quantum con- 
finement, which favors spatial variations of the order parameter, compensating the corresponding 
energy cost. The size-dependent pinning of the chemical potential further increases the spatial vari- 
ation of the pair condensate. The role of nonuniform pairing is smaller in less symmetric confining 
geometries and/or in the presence of disorder. However, it always remains of importance when the 
energy spacing between discrete electron levels 5 is approaching the scale of the bulk gap As, i.e., 
5 > 0.1-0.2 As. 

PACS numbers: 74.20.Fg, 74.78.Na 
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I. INTRODUCTION 

Quantum confinement plays a fundamental role in su- 
perconductors with nanoscale dimensions. Interplay of 
quantum confinement and pairing correlations results 
in important qualitative changes in the superconduc- 
tor characteristics Because of technological reasons 
quasi-OD superconducting structures (i.e., ensembles of 
small grains) were the first where this interplay was inves- 
tigated experimentally. Initial attempts by Giaever and 
Zeller at the end of 60s used tunneling studies on large en- 
sembles of superconducting particlesJ^ Since that time 
most of the studies on superconducting correlations in 
grains were performed with grain powder o^^'^'-' or on films 
made of crystalline granules separated by amorphous 
inter-granular spacCi ^^'^^ In the pioneering work of Ralph 
et a L ^'Si ; .24 the discrete electron spectrum was measured 
for a single grain. Their technique (single-electron tun- 
neling spectroscopy) enabled them for the first time to 
probe superconducting correlations in an individual Al 
grain. Very recently, STM was used to detect the su- 
perconducting gap of an isolated ultra-small lead grain 
deposited onto a sihcon substrate (see e.g. Refs. [25| and 
\2w . These advances opened new prospects to examine 
superconductivity in individual metallic nanograins with 
unprecedented detail, e.g., to investigate the influence of 
the confinement on the superconducting correlations. 

The main feature of a superconducting nanograin that 
makes them different from a bulk superconductor is the 
formation of discrete electron levels with average energy 
spacing 5 « 2Tr'^ h'^ / {mkpV), with kp the bulk Fermi 
wave number and V the system volume. It can be of the 
same order as the bulk gap Ag , or even larger in the case 



of ultra-small nanograins. Therefore, size-quantization of 
the electron spectrum can have a substantial impact on 
the basic superconducting characteristics of such quasi- 
OD superconducting systems. 

The understanding of the fundamental properties of 
superconducting correlations in low-dimensional struc- 
tures, in particular in isolated metallic grains, has experi- 
enced a remarkable development in the last two decades. 
Theoretical aspects, which have attracted the most at- 
tention are the following. The problem of the break- 
down of BCS superconductivity in ultra-small metallic 
grains was addressed in several papersi^i^ii^ The ef- 
fect of the shell structure in the single-electron spectrum 
on the superconducting correlations was pointed out for 
nanograinsSiil and ultrasmall metallic clusters jiSiii The 
ground state properties of the BCS pairing Hamilto- 
nian of ultra-small grains were considered beyond the 
mean-field approximation using the Richardson exact so- 
lutionji"— 

A spatially uniform pairing was assumed in these and 
other works and, as a consequence, the matrix elements 
of the pairing interaction were taken independent of the 
relevant single-electron quantum numbers, i.e., they were 
set to —g/V, with g > the coupling constant and V the 
volumci^ This is, say, a bulk-like approximation recov- 
ered when the single-electron wave functions are taken as 
plane waves. However, the translational invariance is bro- 
ken in nanograins, which leads to a position-dependent 
order parameter. As a result, the pairing gap becomes 
strongly dependent on the relevant quantum numbers, 
which is directly related to a confinement-induced mod- 
ification (as compared to —g/V) of the matrix elements 
controlling the scattering of the time reversed states. An- 
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other important issue is that single-electron levels can 
form bunches and even a kind of shell structure in sym- 
metric confining geometries. In this case the chemical 
potential fi can be pinned to a group of nearly degenerate 
or degenerate levels. This is of importance because the 
density of states in the vicinity of /i strongly influences 
the superconducting correlations. In other words, such a 
pinning plays the role of a filter that selects the contri- 
bution of a particular single-electron shell (or of a group 
of close levels) to the superconducting order parameter. 
Such a contribution is, as a rule, spatially nonuniform. 

The aim of the present paper is to investigate ef- 
fects related to a spatially nonuniform pairing in metal- 
lic nanograins, which was not studied in previous pub- 
lications. For illustrative purposes we consider metal- 
lic spherical nanograins, where the spatial dependence of 
the superconducting condensate is pronounced (the order 
parameter can vary with position by an order of magni- 
tude). In less symmetric confining geometries and/or in 
the presence of disorder spatial variations of the order 
parameter are reduced. However, our study implies that 
nonuniform pairing remains of importance when the in- 
tcrlcvcl spacing S is approaching the scale of the order 
of the bulk gap A b ■ Any remaining grouping of single- 
electron levels, that is always present in real samples, 
even strengthens the effect of interest. We work in the 
mean-field approximation and, thus, stay in the regime 
6 < Ab- Below we consider Sn and Al with Ab ~ 0.616 
and 0.25 meV, respectively (for the parameters used be- 
low). Using the above values of As, we find that the 
mean-field approach is valid for D > 6-8 nm. with D the 
sphere diameter. 

Our paper is organized as follows. In Sec. |lTl we out- 
line the formalism how to obtain a self-consistent solution 
to the problem. In Sec. IIIIl we present our numerical re- 
sults. In particular, in Scc. lIII Al we investigate the effects 
of quantum confinement on pairing correlations through 
the modifications of the matrix elements of the pairing 
interaction and the size-dependent pinning of ^ to single- 
electron shells. Sec. IIIIBI is focused on a spatial distri- 
bution of the pair condensate and its relation to modi- 
fications of the matrix elements and the size-dependent 
pinning of /i. In Sec. IIII CI we discuss the interplay of 
Andreev reflection with quantum confinement, resulting 
in the formation of Andreev-type states and significant 
dependence of the pairing gaps on the relevant quantum 
numbers. A short summary and discussion are given in 
SecHVl 



II. FORMALISM 

The reduction of the system to the nanometer scale 
leads to the formation of a discrete electron spectrum. 
Moreover, in the presence of quantum confinement, the 
translational invariance of the system is broken, and the 
superconducting order parameter is position dependent, 
i.e., A = A(r). For the mean-field treatment of such 



a situation, it is appropriate to use the Bogoliubov-de 
Gennes (BdG) equations^^i^ which can be written as 

E,\u,) = He\u,) +A\v,), (la) 

E,\v,)^A*\u,)-H:\v,), (lb) 

where Ei stands for the Bogoliubov-quasiparticlc (bo- 
golon) energy, A = A(r) (with r the position opera- 
tor) and the single-electron Hamiltonian is referred to 
the chemical potential /i, i.e., 

i?e(r) = ^ + n?)-M. (2) 

We remark that any magnetic effects are beyond the 
scope of the present paper. For simplicity, the confin- 
ing interaction V^(r) is taken as zero inside the specimen 
and infinite outside: V{r) = Vb i^iR — p) with the bar- 
rier potential Vb ^ oo {R = D/2 and p is the radial 
coordinate for the spherical confining geometry). 

As a mean-field approach, the BdG equations should 
be solved in a self-consistent manner 

A(r)=g^(r|^.,)(«,|r)tanh(^), (3) 

i 

where g > is the coupling constant for the effective 
electron-electron interaction approximated by the delta- 
function potential, i.e., (r, r'|$|r, r') = —g5{T — r'). The 
sum in Eq. ([3]) runs over the states with the single- 
electron energy 

= [{u,\He\u^) + {v,\He\Vi)] £ [-h^D, fnOo], (4) 

with ujd the Debye frequency. As is known, the solu- 
tion of the BdG equations has two branches: (i, -I-) and 
(i,— ) (see Ref. l33l ) for which we have Ei^^ > and 
Ei^^ < 0. The sum in Eq. (jS]) should be taken over the 
physical states [the (i, +) branch], i.e., Ei = Ei^+. 

For a given mean electron density the chemical po- 
tential /i is determined from 
2 

ne = y^[fi{ui\ui) + (1 - ft){vi\vi)\, (5) 

i 

with V = l^ri?"^ the volume of the spherical grain. For 
conventional superconductors the energy gap is typically 
much smaller than the chemical potential. As a result, 
/i stays nearly the same when passing from the normal 
state to the superconducting one.— Therefore, one can 
solve Eq. ([5]) in the absence of superconducting order 
(A(r) = 0).^ 

In a spherical nanograin, because of symmetry reasons, 
the order parameter depends only on the radial coordi- 
nate, i.e., A = A{p). Therefore the pscudospinor in the 
particle-hole space can be characterized by the quantum 
numbers of the angular momentum, i.e., {l,m). The an- 
gular part of the pseudospinor '^i is given by the spherical 
harmonics Yim{6^Lp) in polar coordinates {p^6,ip), i.e., 

(ri*,)-y.™(^,^) (:;;[;]), (6) 
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where i = {j,l,m}, with j the radial quantum number 
associated with the quantum-confinement boundary con- 
ditions 



where 



Uji{p)\p=B, = Vji{p)\p=R = 0. 



(7) 

To solve the BdG equations (|la|) and (jlb| numerically, 
Uji{p) and Vji{p) are expanded in the eigenfunctions of 
the single-electron Hamiltonian He [see Eq. ([2])]. In ad- 
dition, iterations should be invoked, to account for the 
self-consistency relation given by Eq. This program 
is significantly simplified by keeping only the pairing of 
the time-reversed states,— which is a standard approxi- 
mation for the problem of superconducting correlations 
in nanograins. In the framework of the BdG equations 
this can be done through the so-called Anderson approx- 
imate solution for which the particle- and hole-like wave 
functions are assumed to be proportional to the single- 
electron wave function. It means that 



Ujiip) = UjiXji{p), Ujiip) =VjiXji{p), 



(8) 



with the radial part of the single-electron wave function 
given by 

= P3/2 ^ / T-^Kaji-s)' (9) 

R-^l^]i+i{aji) R 

with ji{x) the Z-order spherical Bessel function of the 
first kind and aji its j-node. The coefficients lAji and Vji 
(taken as real) obey the standard constraint (see, e.g., 
Refs. [li) 

U^i+Vl = l. (10) 

Then, inserting Eq. ([8]) into Eqs. pa|) and (flbl wc find 
the following set of coupled equations (here Ejim ~ Eji 
and ^jim = 



with 



R 

A,7= J dpp'xUp)^(p) 



(11a) 
(lib) 



(12) 



M. 



'f^ I^PP^ Xyi'ip) X%{p)- 



To derive Ec^. psp . one should keep in mind the prop- 
erty of the spherical harmonics X]ni=-; \^im{Q,^)\'^ ~ 
We remark that Mjniji is nothing else but the 
pairing-interaction matrix element (i', i) (with i = 

{j,l,—m}) averaged over the states with m = —I,... I 
and m' ~ — . . .1' ^ i.e.. 
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(2Z' + 1)(2? + 1) 



E E( 

rn' — — I' 771 — — I 



As seen from Eq. ([T2|) . a spatially uniform order pa- 
rameter means that the pairing gaps A^; do not depend 
on the quantum numbers j and I. This is compatible 
with Eq. (jlSp only when Mjniji docs not depend on j' 
and I' . According to the definition given by Eq. (|T6|) . we 
have Mjniji = Mjijur and, so, if Mj>i>ji does not depend 
on j', I', it does not depend on j, I either. So, we arrive 
at the standard simplified approach of investigating the 
pairing correlations in metallic grains (see the discussion 
in the Introduction). Below we show that the spatial 
dependence of the order parameter can not be ignored 
in superconducting nanograins, which implies significant 
variations of the matrix elements and pairing gaps with 
the relevant quantum numbers. After a numerical solu- 
tion of Eq. (jisp , the position-dependent order parameter 
can be calculated from 



A(p) = ^A(^')(p), 



(16) 



with the shell-dependent contribution A'-''^(p) given by 
A0')(p) = A(2Z + l)4£i^tanh(^). (17) 



and 



2me i?2 



(13) 



A nontrivial physical solution of Eqs. (jllap and (jllbp 
exists only when 



(14) 



The Anderson prescription about the pairing of the 
time-reversed states allows one to rephrase the self- 
consistency relation [see Eq. (jS])] as follows: 

^■^^■'''■'"'^■'"':tanh(^), (15) 
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III. DISCUSSION OF RESULTS 

A. Enhanced intrashell matrix elements and 
quantum-size pinning of the chemical potential 

Numerical calculations were performed with the set 
of parameters typical for ti n'^^i^^ : fiw^i/fc^ = 195 K, 
gN{Q) ^ 0.25, with A^(0) the bulk density of states 
at the Fermi level (we use the bulk electron density 
He = 148nm~3, see, e.g., Ref. l37f). 

Figure [TJa) shows the critical temperature (in units of 
the bulk critical temperature Tc^b) versus the nanograin 
diameter D as calculated from Eq. p5l) when the matrix 
elements of the electron-electron interaction and the size- 
dependent variation of the chemical potential have been 
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FIG. 1. Critical temperature versus tiie grain diameter as 
calculated for: (a) Mjni ji ^ —g/V and /is; (b) Myii ji — 
—g/V and fi / /is; and (c) Mjni^ji / —g/V and fi = /is- 
Tiie dasiied curves in (a) show approximate lower and upper 
boundaries for the quantum-size oscillations of Tc, both curves 
represent the same dependence Tc/Tc,b = 1 + a{D / Dq)^^'^ , 
with Do = 50 nm and a = 1 (the lower boundary) and a — 
3.5 (the upper one). The same curves are also given in (b) 
and (c), for comparison. 



fully taken into account. Results in Fig. [T] are presented 
for a step AR = 0.01 nm. For each radius the critical 
temperature was defined as the temperature above which 
the spatially- averaged order parameter (A(p)) becomes 
smaller than 0.01 of its value at T = 0. Our numerical 
results exhibit two features typical of the size-dependent 



pairing characteristics in high-quality superconducting 
nanograins and nuclei. First, we observe an overall in- 
crease of Tc with decreasing D (it is very pronounced due 
to the highly-symmetric confining geometry). Second, Tc 
oscillates wildly with D. This oscillatory behavior can be 
understood in the following way. The pair correlations 
are nonzero only for the states within a finite range (the 
Debye window) around the chemical potential ^. More- 
over, the main contribution to the sum in Eq. (|T5|) comes 
from the states in the very vicinity of the Fermi level, 

because in this case the expression Aji/ ^S,'^^ + A^j ~ 1 

(^j-; ~ 0). When varying the nanograin size, the num- 
ber of states in the Debye window changes. The smaller 
the diameter, the smaller the number of relevant states 
contributing to the pairing characteristics and, as a re- 
sult, the more significant is such a change. This change 
is not monotonous but rather oscillating due to a per- 
manent competition between incoming and outcoming 
states. As a consequence, all basic pairing characteris- 
tics, e.g., Tc and pairing gaps Aj-;, exhibit quantum-size 
oscillations. It is not only typical of nanograins with su- 
perconducting correlations (see, e.g., the recent paper—) 
but it is also present in superconducting nanowiresi^riS 
and nanofilmsj^i^ Such oscillations are pronounced for 
small diameters/thicknesses but decay with increasing 
the characteristic size so that Tc approaches the bulk crit- 
ical temperature Tc,b (for our parameter T^b = 4.01 K). 
It is interesting to note that the overall increase of Tc with 
decreasing D in Fig. [Tfa) is similar to a size-dependent 
enhancement of the pairing gap in nuclei, where it is pro- 
portional to 1/VA (see, e.g., Ref. [iol ). with A the num- 
ber of nucleons. In particular, the two dashed curves in 
Fig. [ija) show approximate upper and lower boundaries 
for Tc, highlighting the magnitude of the quantum-size 
oscillations: both curves represent the same dependence, 
i.e., Tc/Tc,B = 1 + a{Do/Df/^, with Dq = 50nm and 
a = 1 and 3.5 for the lower and upper boundaries, re- 
spectively [{Do/Df/'^ cx A^e"^^^ with Nc = UeV the 
number of electrons]. We remark that real samples ex- 
hibit inevitable shape and size fluctuations that affect the 
high-degeneracy of single-electron levels. Hence, mea- 
surements on an ensemble of nanograins will significantly 
smooth the quantum-size oscillations in the critical tem- 
perature and reduce its overall enhancement with de- 
creasing nanograin size (see, also. Sec. lIVp. For instance, 
in experimentally fabricated tin nanograins of a semi- 
spherical shape the observed enhancement of the excita- 
tion gap over its bulk value is aboutSI 60% for the particle 
heights « 10-20 nm. This is significantly smaller than the 
enhancement of Tc shown in Fig.[IJa). However, detailed 
investigations of the enhancement of Tc in superconduct- 
ing nanograins is beyond the scope of our present paper. 
Here we arc interested in a spatially nonuniform distribu- 
tion of the pair condensate which is of importance even in 
the presence of shape and size fiuctuations and disorder 
(see the discussion in Sec. IIV[) . 

In order to outline the role of the matrix elements 
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Mjniji [see Eq. p^ ] of the electron-electron interac- 
tion we also show what happens when the true matrix 
elements arc simply replaced by those of the bulk-like 
form: Mjniji = —g/V, which is what is usually done 
when investigating the superconducting correlations in 
nanograins. The results are displayed in Fig.[ljb) and, as 
seen, the difference with respect to Fig.[IJa) is significant. 
To simplify the comparison, we show also in Fig.[IJb) two 
solid curves that represent the radius-dependent upper 
and lower values of T^, from Fig. [IJa). 

TABLE I. Matrix elements Mjniji = Mji^yi, in units of -g/V 
calculated at D = 7.1 nm for quantum numbers such that 
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To clarify the physical reason why using the true ma- 
trix elements leads to significant deviations from the re- 
sults found for Mjniji = —g/V, we show in Table I the 
numerical values of Mjniji (calculated in units of ~g/V) 
for D = 14.2 nm (only the states within the Dcbye win- 
dow are given). As seen, the diagonal (intrashell) ma- 
trix elements Mji^i are strongly enhanced as compared 
to —g/V. However, the matrix elements controlling the 
scattering of the time reversed states between different 
shells (intershell) arc often decreased in absolute value 
with respect to —g/V. So, the question arises why the 
superconducting correlations are enhanced for the true 
matrix elements? The point is that the intershell inter- 
actions are of less importance due to a size-dependent 
pinning of the chemical potential to the groups of de- 
generate or nearly degenerate levels (shells can be often 
close to each other in energy), see the next paragraph. 
When fi is pinned to a particular shell, then the single- 
electron energy measured from fi is zero for the states 
from this shell. These states make a major contribution 
to superconducting correlations unless diameters are not 
large enough D < 20-30 nm, in other words, the num- 
ber of contributing shells is less than 10-15. In this case 
the superconducting correlations are nearly determined 
by the pairing gap Aji associated with the shell pinned 
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FIG. 2. (a) Size variations of the chemical potential, accompa- 
nied by an overall shift of fi to upper values with decreasing D. 
(b) Details of the quantum-size pinning of jj, (filled squares) 
to the single-electron levels (solid curves), small diameters are 
shown for simplicity. 

to fi. From Eq. (|15p it is seen that Aji for the states 
with ^ji = is mainly governed by the intrashell matrix 
element Mjiji. For instance, when ignoring the contri- 
bution of all other states one simply obtains (at T = 0) 

A,7~-(/ + ^)M,,,,,. 

When the diameter increases beyond 20-30 nm, then the 
intershell matrix elements approach —g/V while the in- 
trashell matrix elements are still significantly different 
from the bulk- like behavior. However, the role of the 
states with ^ji = is becoming less and less important for 
larger diameters due to the presence of larger and larger 
number of shells making a contribution to the pairing 
correlations. As a consequence, the difference between 
the data in Figs. [IJa) and (b) decreases when approach- 
ing D — 35-40 nm. together with the amplitude of the 
quantum-size oscillations of Tc. 

In the fully self-consistent scheme the chemical poten- 
tial is determined in such a way that the mean elec- 
tron density rie is constant [see Eq. (O]. However, size- 
dependent variations of ^ are of importance not only be- 
cause they simply prevent the mean electron density from 
deviations. In fact, such deviations are almost insignifi- 
cant: our calculations for ^ = fj,B show that Ue decreases 
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FIG. 3. Spatial distribution of the pair condensate in spherical nanograins: A(p) (calculated at T = 0) versus p for diameters 
= 12nm (a), 13.52 nm (b), 14.2 nm (c), 16nm (d), 16.4nm (e) and 17.54 nm (f). 



by a few percent when D reduces to 10-20 nm. A more 
interesting thing is that the size-dependent variations of 
/i have a pronounced effect on the superconducting cor- 
relations. In particular, this can be seen from Fig. [SJc), 
where Tc is calculated for the true matrix elements and 
fi = fiB- What is the reason for this suppression of Tc? 
In the presence of the formation of strongly degenerate 
electron levels or bunches of electron levels with almost 
negligible spacing between them, the chemical potential 
lies mostly at the highest partly-filled degenerate level 
(see, e.g., Refs. fiol and [TT[). Pairing correlations are sig- 
nificant only within the Debye window around the chem- 
ical potential fj, and are strongest^i exactly at fj,. Hence, 
when fj, is pinned to a shell level, this favors the pair- 
ing correlations at this level and, in turn, through the 
self-consistency relation, favors the pairing correlations 
at neighboring shells. In other words, if the level to which 
the chemical potential is pinned is highly degenerate than 
the phase space for the strongest pair scattering is en- 
larged and, consequently, the system gains in interaction 
energy and, as a result, superconducting correlations are 
strongly enhanced. It is different when fi is not pinned to 
a shell, which is mostly the case for a constant chemical 
potential, e.g., for /i = /i^. Here the relevant shells en- 
tering the Debye window arc as a rule specified by ^ji 7^ 
and, so, their contributions arc diminished. 

The above discussion is further illustrated by our nu- 



merical results for /i in Fig. ^ As seen from panel (a), 
when keeping the electron density of the system con- 
stant, ^ slightly shifts systematically up with decreas- 
ing D and exhibits size-dependent oscillations, as seen 
from Fig. [H^a) . These oscillations are a signature of the 
size-dependent pinning of /j, to groups of degenerate or 
nearly degenerate single-electron levels. This is clearly 
seen from Fig.[2)^b). where variations of /i (filled squares) 

are plotted versus D together with the single-electron cn- 

2 ^2 

ergies measure from the band bottom, i.e., 2^"^ (solid 
curves). For the sake of simple illustration, panel (b) 
shows the data for extremely small diameters, where the 
energy spacing between the shell levels is pronounced 
and, as a result, the size-dependent oscillations of ^ are 
not so wild as it happens for higher diameters. As follows 
from Fig. [IJb) /i is pinned to a shell level in most cases, 
which, as mentioned above, represents incomplete shells. 
Sometimes fi can be found between two neighboring shell 
levels, which corresponds to the case of a fully occupied 
lower shell. 



B. Spatially nonuniform pair condensate 

In the previous paragraph we considered the effect of 
quantum confinement on pairing correlations through the 
matrix elements and quantum-size pinning of fi. As dis- 
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cussed at the end of Sec. [Ill a framework which incorpo- 
rates both issues appears to be only consistent when the 
position-dependent superconducting order parameter is 
taken into consideration. Thus, our resuhs discussed in 
the previous section suggest that the spatial variations 
of A(p) will be pronounced even in nanograins with di- 
ameters up to D = 20-30 nm. However, it is usually 
argued that spatial variations of A(p) cost significant ex- 
tra energy and, so, they are strongly suppressed when 
D with £ the bulk coherence length (see, for in- 

stance, Ref. \2a) . In addition, D should be larger than 
Xp: in practice, kpD ~ 10 is assumed to be sufficient 
to ignore any spatial dependence of the order parame- 
terpi^iii For typical metallic parameters kpD ^ 200-400 
for D = 10-20 nm and this is the reason why the spatial 
dependence of the order parameter was ignored in most 
papers on superconducting correlations in nanograins. 
To go in a more detail on this point, we below discuss 
our numerical results on A.{p). 

In Fig. [3] the radial dependence of the superconduct- 
ing order parameter is shown as calculated from Eq. ()16|) 
for D = 12 nm (a), 13.52 nm (b), 14.2 nm (c), 16 nm (d), 
16.4 nm (e) and 17.54 nm (f). The shells making a contri- 
bution to the superconducting correlations are also dis- 
played in each panel, and the quantum numbers of the 
shell level pinned to fj, are underlined. As seen, we in 
general have a nonuniform distribution of the pair con- 
densate for diameters D ~ 10-20 nm, which is in agree- 
ment with our expectations. For example, let us consider 
the results plotted in panel (c). Here ^ is pinned to the 
shell level = (101,1) and, so, single-electron states 
with J = 1 and I = 101 make a major contribution to 
A(p), which results in a significant enhancement of the 
order parameter next to the edge, i.e., for p/R = 0.9-1.0. 
The profile of this enhancement is determined by the ra- 
dial wave function Xi loi(p) with two pronounced local 
maxima (A/As = 14'.3 and 7.2 at p/R = 0.9 and 0.97, 
respectively) and one node (recall that j is the number 
of the nodes of the radial wave function). All the other 
shells displayed in Fig. ^c) are specified by ^ and, 
as a result, their contributions is much less significant. 
The local maximum A{p)/Ab = 2.3 at p/R^ 0.1 is due 
to states (j, I) = (31, 11). The shells with (j, I) = (23, 29) 
and (19, 29) are responsible for local enhancements of the 
order parameter up to 1.5-2. OAs at p/R ~ 0.27 and 0.36, 
respectively. At last, the shell (8, 71) produces the local 
maximum at p/R = 0.64. In general, the larger the an- 
gular momentum, the larger the values of p/R at which 
the corresponding single-electron states have an effect on 
the profile of A(p). 

It is worth noting that typically, the order parame- 
ter is strongly suppressed in the center (p = 0) except 
of rare cases when states with zero angular momentum 
contribute to the pairing correlations. One such exam- 
ple is given in Fig. [^d), where a narrow pick can be 
seen at p = due to the contribution of the shell with 

ao = (4i,o). 

From Fig. [3] it follows that the radial distribution of 



the pair condensate remains strongly nonuniform even 
for D 20 nm. We would like to note that when select- 
ing concrete values of D for Fig. [31 we did not even take 
diameters for which Tc is close to the upper dashed curve 
in Fig. [Tfa). In the case of a strong enhancement of Tc 
the radial distribution of the pair condensate is as a rule 
strongly nonuniform. The points selected for Fig. |3] are 
mainly in a vicinity of the lower dashed curve in panel (a) 
of Fig. [TJ for D = 13.52, 14.2, 16 and 16.4nm we have 
Tc/Tc,B = 6.41, 6.48, 4.082 and 3.78, respectively. How- 
ever, even in this case the order parameter can vary with 
position by an order of magnitude. Spatial variations of 
A(p) are significantly relaxed only when D approaches 
30-40 nm, as seen from Fig. [H 

For our parameters kp — 16.4 nm and, so, we ob- 
tain kpD « 300 for D « 20 nm. Hence, the criterion 
kpD 3> 1 is not very useful in order to estimate the ef- 
fect of spatial variations of the pair condensate. Based 
on our numerical study, we would like to suggest another 
criterion related to a more sensitive energy scale, which in 
the superconducting state is governed by the bulk pairing 
gap A B ■ The spatial distribution of the order parameter 
is always strongly inhomogeneous when S ^ Ab (here it 
is even better to replace As by the size-dependent pair- 
ing gap). The spatial variations decay with a decrease in 
the ratio of the mean interlevel spacing to the bulk order 
parameter, i.e., S/Ab, and our numerical results sug- 
gest that such variations are significantly reduced only 
when S/Ab < 0.05-0.1 (recall that effects of a magnetic 
field are beyond the scope of our paper) . For Sn spheri- 
cal superconducting grains this regime is achieved when 
D > 40-50 nm (note that S sa 2T:'^h^ / {mkpV) underesti- 
mates the intcrshell spacing for spherical confining poten- 
tial) . Despite that our results are for a highly symmetric 
confining geometry, we can expect that the order param- 
eter will be always spatially nonuniform for 5/ Ab > 0.1, 
even when shape imperfections and disorder dissolve a 
shell structure. The reason is that the number of con- 
tributing states (i.e., the states in the energy interval 
~ [p, — Ab,P + Ab) is not very large for S/Ab > 0.1. 
In this case the states pinned to p always make a major 
contribution to the order parameter and, so, the profile 
of the squared absolute value of the corresponding wave 
function will mainly determine the spatial distribution 
of the condensate. Thus, the domain S/Ab ~ 0.1-1.0 
is in general characterized by strong effects due to the 
spatially nonuniform pairing. 

We remark that our conclusions do not contradict the 
usual argument that spatial variations of the order pa- 
rameter cost extra energy. Let us compare a bulk super- 
conductor with a superconducting nanograin. In bulk the 
relevant matrix elements controlling the scattering of the 
time- reversed states are —g/V and the order parameter is 
spatially uniform (in the absence of a magnetic field). As 
opposed to bulk, the pair condensate significantly varies 
with position in nanograins, which results, of course, in 
an increase of the kinetic energy. However, the intrashcll 
matrix elements are now enhanced in absolute value as 
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FIG. 4. The order parameter A(p) for sufficiently large diameters D — 33.6 nm (a), D — 35 nm (b) and D = 36 nm (c) 



compared to —g/V due to quantum confinement. This 
compensates energy costs of spatial variations of the or- 
der parameter. 

The discussion in the previous paragraph is also related 
to arguments that invoke the conventional Ginzburg- 
Landau theory. According to this arguments the or- 
der parameter is uniform in samples with size smaller 
than the bulk coherence length. When applying this 
to nanograins, one can conclude that the pair con- 
densate should not vary with its position. However, 
this is not true. It is well-known that one should 
be careful when applying the conventional Ginzburg- 
Landau theory to superconductors with characteristic 
size smaller than the zero-temperature (BCS) coherence 
length Strictly speaking, Gor'kov's derivation of 

the conventional Ginzburg-Landau formalism from the 
BCS approach is not applicable on a scale smaller than 
^0 (see,.e.g^ Ref. HI). For Sn we have £_q « 230 nm (see, 
e.g., Ref. I32I ). Thus, in the case of interest I? <C Co, and 
one can hardly invoke the conventional Ginzburg-Landau 
formalism to check whether or not A(p) varies with p. 



C. Confinement-induced Andreev-type states 

Here we would like to discuss one more issue related to 
a spatially nonuniform pairing in nanograins. This is the 
formation of Andreev-type states induced by quantum 
confinemen t ^^'^^ (see also a similar paper— discussing 
Andreev-type states in an ultracold trapped superfluid 
Fermi gas). Since the 60s (see Refs. I42l - l44l ) it is known 
that quasiparticles can "feel" a spatial variation of the 
superconducting order parameter as a kind of potential 
barrier. This physical mechanism (referred to as Andreev 
mechanism below) is the basis for Andreev quantization 
investigated previously for the core of a single vortex for 
the mixed state of a type-II superconductor— and for 
an isolated normal region of the intermediate state of 
a type-I superconductor— (or for a similar case of SNS 
contacts^). Based on our consideration of Sec. lIII Bl one 
can expect that Andreev-type states can play a remark- 



able role in superconducting nanograins due to significant 
spatial variations of the superconducting order parame- 
ter. This is very similar to recently investigated Andreev- 
type states in superconducting nanowires /nanofilm a^^'^^ , 
where the pair condensate is position dependent in the 
direction perpendicular to the nanowire/nanofilm due to 
the quantization of the perpendicular electron motion. 
In Ref. 45 it was shown that 



A,; 



dVA(r) \u,ir)f + \v,{r)\' 



(18) 



which means that the pairing energy gap A^ is the 
averaged value of the order parameter "watched" by 
the quasiparticles with quantum numbers i. Note that 
|ui(r)p -I- |t;i(r)p can be interpreted as the spatial dis- 
tribution of quasiparticles according to the well-known 
constraint / d^r(|ui(r)p -I- |z;i(r)p) = 1 [see, e.g., Ref.lssI 
and Eq. ([TO])]. When inserting Eqs. ^ into Eq. (HHI, 
one can easily obtain Eq. ((T2|) with A^ = Aji. If quasi- 
particles avoid the domains of enhanced pair conden- 
sate, the corresponding integral in the right-hand-side of 
Eq. becomes smaller and, hence, such quasiparticles 
have smaller pairing gaps Aji. They can be referred to 
as Andreev-type states. 

Our numerical study of quantum-number dependent 
pairing gaps Aj-; for metallic nanograins reveals a sig- 
nificant role of Andreev mechanism. Let us consider 
D = 13.52 nm, the corresponding spatial distribution of 
the pair condensate is given in Fig. [3l[b) . To show how 
different species of quasiparticles are distributed in the 
radial direction in this case, the radial-dependent shell 
contributions (at T 0) A(j')(p) [see Eqs. ^ and 
([T7])] are plotted in Fig. EJa). We remark that such a 
representation is more informative than simply a plot 
of \uji{p)\'^ + \vji{p)\'^. First, the radial dependence of 
oc Xjiip) is the same as that of \uji{p)\'^ + \vji{p)\'^ cx 
Xjiip) [see Eq. Second, a plot of A(^')(p) gives also 
information how the corresponding states contribute to 
A(p). From Fig. [3Jc) we can see that a significant en- 
hancement of the order parameter occurs at p/R = 0.45- 
0.7. From Fig. O^a) it is clear that this enhancement 
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is due to the states with {j,l) = (14,48) and (9,63). 
Other shells, i.e., (27, 16) and (23, 25) contribute less, and 
the corresponding quasiparticles, representing Andreev- 
type states, are mainly located beyond the domain p = 

0. 45 — 1.0. As a result, they have smaller pairing gaps, 

1. e., A27a6 = 2.65 As and A23.25 = 2.81 As, as com- 
pared to Ai4,48 = 4.098 Ab and Ag^ea = 5.77 A^. As 
seen, the quasiparticles with (j, /) = (27, 16) are most 
successful in avoiding the local enhancement of A(p) at 
p/R = 0.45-0.7 and, so, A27^i6 is the smallest pairing 
gap. Such a manifestation of Andreev mechanism is not 
a particular feature of = 13.52 nm. In general, Aji 
strongly varies with j and I for diameters < 30 — 40 nm, 
i.e., where spatial variations of the order parameter are 
still pronounced. Quite often such variations can be an 
order of magnitude, as, e.g., for D = 14.2 nm (see A(p) 
given in Fig.[2jc)). At this diameter a great enhancement 
of A{p) takes place at p/R = 0.9. This is due to the con- 
tribution of the shell with (j,/) = (1,101) [see Fig. [Sib)]. 
Other shells make much less important inputs and the 
corresponding quasiparticles are mainly distributed be- 
yond the domain p/R ~ 0.9-1.0. So, as compared 
to Ai^ioi ~ 9.32 As, they have significantly smaller 
pairing gaps, i.e., Asi^n = 1.6Ab, A23,29 = 1.62 As, 
Ai9,39 = 1.72 As and Ag,?! = 2.35 A^. Thus, the inter- 
play of Andreev mechanism and quantum confinement is 
responsible for variations of Aji with the relevant quan- 
tum numbers. 

One could expect that such a serious difference in pair- 
ing gaps of different quasiparticle species can result in a 
pronounced drop of the ratio of A^; (the minimal en- 
ergy gap) to the critical temperature ksTc, similar to 
the case for quantum superconducting nanowiresJ^ The 
main idea here is that A^; is governed by Andreev- type 
states and, hence, is decreased. Unlike A^;, Tc is con- 
trolled by the quasiparticles making a major contribu- 
tion to A.(p) and, so, Tc is coupled to their higher pair- 
ing gaps. As a result, AE/ksTc can be significantly 
smaller than in bulk. For instance, one can expect 
that Ae = Aai^ii = 1.6 Ab at D = 14.2 nm while 
Tc is governed by Ai^ioi ~ 9.32 A^. However, this is 
not correct for nanograins. The point is that A^; is a 
spectroscopical gap which is probed by STM. It is de- 
fined as Ae = mhijiEji. For nanowires the subband- 
dependent pairing gap is always the minimal quasiparti- 
cle energy due to a quasi-free spectrum in the direction 
parallel to the nanowirc. For nanograins this is differ- 
ent. In particular, for D = 14.2 nm we have the fol- 
lowing single-electron energies (absorbing p) of the rel- 
evant shells: ^31.11 = -18.6 Ab, 63,29 = -26.03 As, 
^9,39 = 20.9 As, 6,71 = 22.6 As and 6,101 = 0. Hence, 
one can calculate that A^; = i^i.ioi = Ai,ioi in spite 
of the fact that Ai,ioi is the largest pairing gap. Thus, 
although Andreev mechanism plays a significant role in 
superconducting nanograins, it can hardly be probed by 
STM-mcasurcmcnts due to the nonzero interlevel spac- 
ing, unlike quantum superconducting nanowires. 



Z)=13.52 nm 

(/J)=(27,16) (9,63) 
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FIG. 5. Shell-dependent contributions to the order param- 
eter Aji(p) for relevant shells: (a), D = 13.52 nm, {j,l) — 
(27, 16), (23,25), (14,48) and (9,63); (b) D = 14.2nm, 
(i,0 = (31,11), (23,29), (19,39), (8,71) and (1,101). 



IV. CONCLUSIONS AND DISCUSSION 

In conclusion, we have shown that the spatial distribu- 
tion of the pair condensate is essentially nonuniform in 
metallic nanograins. In particular, the spatially nonuni- 
form pairing can proliferate in nanograins even when 
kpD ~ 300 and, so, the usual criterion to neglect vari- 
ations of the superconducting condensate with position, 
i.e., kpD ^ 1, is not very useful and can result in wrong 
conclusions. This is the reason why effects due to spa- 
tially nonuniform pairing in superconducting grains were 
previously overlooked. Our study suggests that a new 
criterion should be based on a more delicate energy scale 
(as compared to the Fermi energy), which, in the su- 
perconducting state, is given by the bulk order parame- 
ter Ab- It turns out that the pairing becomes spatially 
nonuniform when the interlevel spacing S exceeds 0.1- 
0.2 As. Variations of the order parameter with position 
exhibit a pronounced enhancement with an increase of 
5/Ab- When 5 ^ As, such variations can be almost 
an order of magnitude in highly symmetric grains. At 
first sight, this seems impossible because it costs extra 
energy for such spatial variations. However, a nonuni- 
form distribution of the pair condensate is accompanied 



10 



by enhanced pairing interaction matrix elements, which 
compensates the energy cost for an inhomogeneous dis- 
tribution of the condensate. Another point is the size- 
dependent pinning of the chemical potential to groups 
of degenerate or nearly degenerate energy levels. Such a 
pinning plays the role of a filter that increases the contri- 
bution of the single-electron levels in the vicinity of the 
chemical potential and suppresses contributions of other 
states. This results in an additional mechanism favoring 
spatially nonuniform pairing in metallic nanograins. 

In this paper we investigated a highly symmetric con- 
fining geometry. Due to this feature the problem becomes 
effectively one-dimensional (the order parameter depends 
only on the radial coordinate) and, so, sufficiently large 
diameters up to D ~ 40 nm can be investigated. This 
size is almost impossible to reach theoretically for grains 
with the order parameter depending on three relevant co- 
ordinates due to time consuming numerical calculations. 
Such an effectively one-dimensional problem has large de- 
generation factors for the corresponding shell structure, 
resulting in a significant enhancement of the pairing cor- 
relations. In reality there can be several issues that may 
lead to a splitting of the shell levels. It will decrease the 
degeneration factors and, so, reduce the pairing correla- 
tions, since the main contribution to the sum in the gap 
equation comes from the transitions within the same shell 
pinned to the chemical potential. Among such issues is 
the Jahn- Teller deformation, i.e., the transformation of 
a spherical nanograin with incompletely filled shells to 
an ellipsoidal shape. In addition, the surface imperfec- 
tions and impurities can significantly change the distri- 
bution of single-electron levels. However, our qualitative 
results arc quite generic and do not depend on a par- 
ticular shape of a nanograin and the presence of possi- 
ble imperfections. For instance, when i5 '--^ the pair 
condensate will always be spatially nonuniform because 
only a few single-electron levels enters the energy interval 
« — Ab, ^ + Ab]- Due to the dominant contribution of 
such levels to A(r), one can expect that the pair conden- 
sate acquires a profile governed by the squared absolute 
value of the wave function for the single-electron state 
closest in energy to /i. This is significantly strengthened 
by an increase (in absolute value) of the diagonal matrix 
elements and, in addition, by the pinning of 

the chemical potential to the single-particle levels. 

We remark that the diagonal matrix elements, i.e., 
[see the definition for $ below Eq. ([3])] are al- 
ways enhanced as compared to —g/V in the presence of 
quantum confinement, whatever disorder and shape im- 
perfections. This can be seen from the following simple 
arguments. Introducing </5i(r), the wave function associ- 
ated with state i, one can write 
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FIG. 6. (a) Single-electron energies (given in units of 
the Debye energy hwo) ordered in ascending manner versus 
the ordering number for the rectangular-shaped aluminum 
nanograin with Lx ~ 7.06 nm, Ly — Lx/1-1, = 1.1 Lx 
(squares) and for a cubic nanograin with Lx ~ 7.06 nm 
(triangles), (b) The size-dependent excitation gap Ab/As 
versus Lx (in steps of SLx ~ 0.01 nm) for an ahiminum 
nanograin of the rectangular shape with the dimensions 
Lx,Ly = Lx/l-l,Lz = 1.1 Lx'- squares represent the re- 
sults calculated with the modified matrix elements and with 
proper variations of ji; stars are the data obtained for the 
bulk-like matrix elements —g/V and jj, = /is. (c) The spa- 
tial distribution of the pair condensate in the rectangular 
grain with Lx — 7.06 nm and Ly — Lx/1-1, Lz = 1.1 Lx, 
A[x) = I^{x,y,z)\y^x/i.i,z=i.ix- 



element can be rearranged as 



1 + V d^rdfir) 



Due to the normalization condition we have |(pi(r)| = 
y + di{r), where / d'^rdi{r) — 0. Then, the above matrix 



The second term in the brackets is always positive in the 
presence of quantum confinement, i.e., when di(r) ^ 0. 
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It is zero only when (/5i(r)'s are chosen in the form of 
plane waves, which results in = —g/V. 

The above discussion can be supplemented by our nu- 
merical results calculated from the BCS-like equation 
similar to Eq. (jlSp but now for aluminum nanograins 
of rectangular shape with dimensions L^^ Ly = 
L^/l.l, = I.IL^,. For aluminum we have^i^i^: 
HuJo/kB = 375 K, gN{0) = 0.18, and hb = 
11.67cV, which corresponds to the electron density 
Ue = 181 nm"-^. In Fig. [U^a) single-electron levels ar- 
ranged in the ascending order arc shown within the De- 
bye window for the rectangular nanograin with Lx ~ 
7.06 nm (squares). The same is also given here for a cubic 
aluminum nanograin with = Ly = Lz = 7.06 nm (tri- 
angles). As seen, single-electron levels for the rectangu- 
lar shape arc distributed in a nearly equidistant manner 
(with S w 0.2-0.3 mcV ^ Ab) contrary to the states cor- 
responding to the cubic geometry. It is well-known that 
an almost equidistant distribution^ of single-electron lev- 
els near ^ is also expected in the presence of significant 
imperfections such as the surface roughness and/or im- 
purities. So, our results in Fig. |6] give a feeling about 
the role of the spatially nonuniform pairing in disordered 
metallic grains. The excitation energy gap A^; for the 
rectangular nanograin is shown in units of in Fig. [6] 
as a function of in the interval = 7-8 nm. Here 
squares represent our results calculated with the modi- 
fication of the matrix elements and with /x varying with 



Lx] stars are the results found for the bulk-like matrix 
elements —g/V and ^ = fiB- As seen, A^ (oc Tc) is 
now two-times enhanced as compared to As (on aver- 
age), which is much less significant than for highly sym- 
metric grains (compare with Fig. [l]) due to a splitting of 
the shell levels. However, the effect of interest is still pro- 
nounced: Ae calculated for the modified matrix elements 
and with account of size variations of ^ is generally larger 
by a factor of 1.5-2.0. The spatial profile of the order pa- 
rameter is nonuniform with local enhancements over its 
average value by about 100%, see, e.g.. Fig. M,c)- For 
rectangular grains with = 7-8 nm we have S ^ Ab- 
However, as we checked, the spatially nonuniform pairing 
and the related effects of the modification of the relevant 
matrix elements and the size-dependent pinning of ^ are 
of significance even for smaller 5's, i.e., when S > 0.1- 
0.2 As {Lx < 14-15 nm). For instance, at Lx — 11 nm 
the order parameter exhibits variations of about 30-40% 
of its averaged value. These results are in agreement with 
our expectations based on the investigation of the highly 
symmetric spherical grains. 
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^® Small corrections to — p- of about 1% — 3% were considered 
in Ref. [l3 for chaotic grains. 

Impact of a shell structure on superconducting correlations 
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